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Abstract 

The modeling and investigation of the dynamics and configurations of animal groups is a subject 
of growing attention. In this paper, we present a continuum model of flocking and use it to 
investigate the reaction of a flock to an obstacle or an attacking predator. We show that the 
flock response is in the form of density disturbances that resemble Mach cones whose configuration 
is determined by the anisotropic propagation of waves through the flock. We analytically and 
numerically test relations that predict the Mach wedge angles, disturbance heights, and wake 
widths. We find that these expressions are insensitive to many of the parameters of the model. 
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I. INTRODUCTION 



The flocking of biological organisms into groups has been a phenomenon of long standing 
interest. Birds, fish, bacteria, and certain robots exhibit rich collective behavior. Research 
in this area has generally employed two main modeling paradigms: discrete individuals and 



continuous densities of individuals 
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While both approaches are useful and contribute to a more complete description of flocking, 
we focus on a continuum approach. In this paper, we will consider the response of a flock 
to a stationary or moving 'obstacle'. In the case of a moving obstacle, our considerations 



might also be considered as modeling the avoidance response of a flock to a predator. 
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research investigating obstacle avoidance has employed a discrete approach 
As compared to a discrete description, the continuum description has the advantage of 
economically treating very large numbers of individuals and is, in some cases, easier to treat 
analytically and to interpret. Its disadvantage is primarily that taking the continuum limit 
is an abstraction from the real case of discrete flock members. 

The current paper introduces a moving obstacle into a large flock and studies the effect 
of this obstacle on the flow around the obstacle. We model the obstacle as a localized region 
exerting a repulsive 'pseudo-force' on the flock continuum. Using our description, we are 
able to describe the propagation of information in terms of a few parameters in the model. 
To do this, we use a fluid characterization of a flock. For a review of this type of approach, 
as well as other approaches to modelling flocks, see 21 1. 



To facilitate our analysis we will utilize a linearized theory in which the flock response 
to the obstacle/predator pseudo-force is assumed to be proportional to the pseudo-force 
strength. That is, the obstacle/predator is treated as a linear perturbation. Results obtained 
through this type of analysis are expected to yield qualitative insights to the dynamics of 
the full nonlinear problem, and may also yield quantitative understanding in the region 
far enough from the obstacle/predator where the perturbations become small. In the next 
section, we will introduce our continuum description of the flock. In the following sections, 
we explore the small amplitude wave dispersion relation and derive an expression for the 
disturbances that propagate through the flock. Next, a linearized response is investigated 
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and the resulting density perturbation is analyzed analytically. Finally, results of numerical 
evaluation of the density perturbation are presented and compared to the theory. 



II. CONTINUUM FLOCKING EQUATIONS 



The equations we consider for flocking in three dimensions are 
dv 1 / v'^\ 1 

^il + vVv = - l-^ v - -VP(p) -VU- W(v) 



dt 
dp 

at 



p 



V ■ (pv) = 0, 



(1) 
(2) 



where p is the number density of flocking individuals, v is the macroscopic vector velocity 
field of the flock, v is its magnitude, and Vq and r are constants. The basic structure of 
these partial-differential equations includes terms that define the acceleration of the fluid 
density of the flock, along with continuity of flock members. The right-hand side of Eq. ([T]) 
consists of four 'pseudo-forces' representing speed regulation, pressure, pairwise attraction, 
and a 'non-local viscosity' term. These terms are discussed below. 

The first term on the right-hand side of Eq. ([T]) acts as a speed-regulation term used 



commonly in the literature 
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22] and apparently first used by Rayleigh 



23| as cited by 



|12l | . This term either increases or reduces the magnitude of the velocity depending on how 
the velocity compares to Vq. If f > Vq, the acceleration is negative in the direction of v, and 
thus |v| = f is reduced. If f < Vq, the acceleration is positive in the direction of v, hence v 
is increased. Thus Vq can be regarded as modeling the average preferred natural speed of an 
individual. The time scale for this velocity clamping is r. Note that this speed-regulation 
term is frame dependent and applies when considering the frame in which the medium (e.g., 
air for birds, water for fish, or land for ungulates), through which the flock individuals move, 
is stationary. 

In order to model the presumed tendency of nearby flock members to repel each other 
to avoid collision, some past models have introduced a pressure-like term, as in the second 
term on the right-hand side of Eq. ([T]). Examples can be found in 2l|. In addition, another 



means to model repulsion is via a general repulsive potential; i.e., a pairwise non-local soft- 



core potential (see 



11|). We model repulsion using a pressure term, P{p). For future 



reference, we write the pressure as a Taylor series around a density po as 



Pip) = ci6p + 



6p' + ... (3) 

P=PO 



dp^ 

where = |^|^^^^ and Sp = p- po- 

The third term on the right-hand side of Eq. ([T]) is a long-range attractive pseudo-force 
where long range attraction is used to model the tendency for flocks to form. This force is 
taken to be due to an attractive pseudo-potential, U, which is of the form 

Ui.)^lui.-.W. (4) 

It proves convenient to choose the kernel ^(x— x') to satisfy the modified Helmholtz equation, 

(V^ - kJ) m(x - x') = 47ruo6 (x - x') , (5) 

where Mq > is the strength of the potential. In three dimensions ^(x — x') has the form of 
an attractive exponentially-screened Coulomb potential, 

m(x - x') = -Mo- — . (6) 

X — x' 



The quantity k ^ provides a long-distance cutoff to the attractive pseudo-force. 



his type 



of attractive potential has been used in previous continuum flocking models 

Similar to the non-local attractive potential, we model the presumed tendency for nearby 
flock members to attempt to align their velocities by use of the term 

W(x) = y u;(x - x') [v(x') - v(x)]dx', (7) 

with the kernel w(x — x') satisfying an equation similar to that for the attractive kernel, 

(V^ - 4) w(x - x') = AttwoS (x - x') , (8) 

with strength wo > and screening length scale k,~^. This term reorients the velocity vector, 
v(x), toward the average velocity of the other flock members, weighting velocities of flock 
members closer to x more strongly than those farther away. These are our general equations 
that model flocking. Various dynamical behaviors and flocking equilibria can be explored 
using this framework. However, in the rest of the paper, we consider perturbations around 
a specific equilibrium density defined below. 



We consider the following simplified situation. A particular, spatially-homogeneous 
steady-state solution to Eqs. ([1]) and ([2]) is 



Alternatively, we may think of this equilibrium as a localized approximation of a more 



middle of a nonuniform flock, the density in equilibrium will be nearly constant (see [5i]). 

To the general equations Eqs. ([T]) and ([2]), we will add an additional, external, localized, 
repulsive potential that we view as modeling the effect of a stationary obstacle or a predator 
moving through the flock with velocity Vp. In the next section, we treat this problem within 
the framework of linearized theory and consider how perturbations propagate through the 
flock. 

III. HEURISTIC DISCUSSION OF INTERACTION WITH AN OBSTACLE 
A. Dispersion Relation and Plane Waves 

By looking at the dispersion relation of linear waves in the full system, we may determine 
how such waves propagate within the flock. This will inform us as to the relationship 
between the frequency, wavelength, and propagation direction of the waves. In our case, we 
will flnd that this will predict a disturbance cone when an obstacle or predator is encountered 
by a flock. For convenience we make a transformation of independent variables such that 
x' = X — \pt and t' = t where Vp is the velocity of the obstacle or predator relative to the 
stationary frame in which the preferred flock speed is vq. This is similar to a Galilean frame 
transformation except that the velocity v remains the velocity in the stationary frame. After 
making this transformation we drop the primes on x' and t'. Writing Eqs. ([T]) and ([2]) in 
this new frame gives 



p(x) = po = const, and v(x) = vq = const. 



(9) 





1 / 

V = - 1 - — 



) 



V 



p 



1 



VP(p) -vu-w 



(10) 



- Vp ■ Vp + V ■ (pv) = 0. 



(11) 
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Setting p = po + Sp and v = vq + Sv, with \Sp\ ^ po and |(5v| ^ vq, we substitute these 
into Eqs. ffTOl) and ffTTl) and only keep hnear terms in 6p and (5v. Taking Fourier transforms 
in both space and time we obtain 

—iuj5v + ik • Vr 5v = ^ (vo • 6v) vq — ik | c^— + tt^p ) — i^wi^v (12) 

^^0 V Po / 

= —iujdp + ik ■ v^, 5p + po^k ■ 5v, (13) 

where we define = vq — Vp, and / = /(k, uj) denotes the Fourier transform of /(x, t) given 
by 

/ = /(k, to) = J /(x, t) exp(zci;t — zk ■ xjdxdt. (14) 



Hence we have 



^(^^) - pf^, (16) 
uUk') = w{e) - w{0) = (^) (17) 



BP 

^ P=PO 

Using Eq. (fT3|) to ehminate 6p/ po from Eq. (fT2|) . we arrive at 

r/ , N ■ T r~ 2ivo • 5v (c^ + poM)k ■ (5v, 
a; - k ■ V, + zz/^ (5v = ^^vq + — k. 19 

Restricting our attention to the case where r — and z/^ ^ (the hmit in which speed 
regulation occurs instantaneously and the non-local viscosity is absent), we obtain partic- 
ularly simple results describing the propagation of waves within the fiock. Note that to 
accommodate the r ^ limit in Eq. f|T9l) . we must have that vq • 5v ^ 0. Without loss of 
generality, we choose vq to be in the x direction, which means that 5\ = 6vy y + 6vz z. If 
we now project Eq. f|T9l) onto the y and z directions, we get two coupled equations for 6vy 
and 6vz. These yield the dispersion relation, 

{LU-k- v,)2 = klcl, (20) 

where we have taken z/^ — 0, and defined k_L = ky-y + k^z giving the magnitude as k\ = 
ky + kl- Also, we have replaced + pou ^ c^, which is true for large k. We can write the 
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final dispersion relation as 

UJ = h ■ Yr ^ k±Cs. (21) 

Thus, the group velocity of waves, in the frame moving at a velocity Vp, within the flock is 
given by 

duj k I , ^ 

In the next section, we use this result to derive a disturbance cone that propagates through 
the flock when the flock encounters an obstacle or predator. 

B. Mach Cones 



Following Mach's well-known construction (see for example [2J]) for the cone produced in 
supersonic velocities through a fluid, we may develop a prediction for the information cone 
that is propagating though the flock using the dispersion relation derived above. In the case 
of a stationary object (Fig. [T](a)), the only way that information can travel is perpendicular 
to the direction of motion with the propagation speed of Cg, in the frame of the obstacle, as 
can be seen in Eq. (122|) with v,, = vq (or equivalently Vp = 0). Accordingly, we get a right- 
circular cone in three dimensions (a wedge in two dimensions) of cone angle 6, measured 
from the direction of the flock, given by 

tan^ = — , (23) 

where is defined above (Eq. f|T8l) ). Notice that this is valid for all velocities, contrary to 
a usual acoustic Mach cone which only exists for velocities of the moving object that are 
above the sound speed. Equation [23] limits to an angle of 6' = | for small Vq, and 6 = for 
large vq Cg. 

For the case of a moving obstacle or predator, Eq. fl22l) implies the construction shown 
in Fig. [D^b). From this construction, we obtain the wedge angles, in a plane passing 
through the cone's axis (defined by vq and Vp), given by 

tan^i = ±^-T ( — ) + tan(^). (24) 



Cross sections for the wedge shapes of both the static and moving obstacle are shown in Fig. 
[TJ In three dimensions, a moving obstacle produces an oblique circular cone as illustrated 



in Fig. [21 If vq and Vp are co-linear, the the cone is a right-circular cone with 6'+ = and 

tlj = 0. 



IV. LINEARIZED THEORY IN A TWO DIMENSIONAL FLOCK 

In order to assess the extent to which these predictions apply more generally, we consider 
a two-dimensional case for both the static and moving obstacle situations. We first solve 
a linearized version of Eqs. (ITOl) and (ITT]) for the density fluctuation Sp/po in terms of an 
integral. We then specialize to a static case to analyze Sp/po analytically. Following that, we 
numerically evaluate our integral-expression solution for both the static and moving cases 




FIG. 1: (color online) Diagram for a density disturbance caused by flock moving past an obstacle in 
steady state, a) Static obstacle: Obstacle is stationary (at point O), the dashed red lines indicate 
the intersection of a plane passing through the axis of the cone of disturbance and a large-density 
fluctuation wake that would exist downstream of the flock/obstacle interaction. The angle is 
referred to as the 'wedge' angle, b) Moving obstacle: In the frame of the obstacle at point O, the 
density disturbance is indicated by the dashed lines, vq is the velocity of the flock and Vp is the 
velocity of the obstacle, both in the frame of the medium. 
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and compare the results to our simple predictions above. 

To specialize to two dimensions, we consider a flock equilibrium that is spatially uniform 
in the z direction. Suppose that an obstacle (call it a predator) is moving through the flock 
at constant velocity, Vp, relative to the flxed frame of the medium (e.g. air or water) in 
which the flock moves. We assume that there is no motion in the z direction. We model the 
obstacle by a moving, localized, repulsive potential, //(x, t), and add the term — Vr] to the 
right-hand side of Eq. (fTOj) . We take 77 to be Gaussian in space and given by 

[x - Vpjy + {y- Vpyty 



r] = r]{x,t) = ?7oexp 



(25) 



where rjo is the strength of the obstacle, / is the length scale over which the obstacle acts, 
and Vpx and the components of is the predator's velocity, Vp, in the x-y plane. In 

two dimensions, this can roughly be thought of as a kind of moving flagpole around which 




FIG. 2: An oblique circular cone 
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the flock must navigate. Without loss of generahty, we set vq = f o x. Given this situation 
we consider the steady-state flock response in the approximation of inflnite flock size. The 
dynamics of a flnite-size flock as it impinges on an obstacle hitting a flock has not been 
considered in the present work. For a simulation of such a situation, see [18,] . 

We add the obstacle potential to Eq. (ITOl) and, in the frame of the obstacle, linearize 
around the constant density, as we did in Sec. 1111 A[ Taking a spatial Fourier transform 
of Eqs. (fTOl) and (fTTi) (including the obstacle), we obtain the following steady-state (i.e., 
d/dt = 0) equations, 

(-ik -Vp + ik- Vq) (5v = — ^ (vq ■ 6v) vq - ik \ c^— + u6p + ) - z/^^v (26) 

\ Po J 

—ik ■ Vp 5p + ik • Vq 5p + pQik • (5v = 0, (27) 

with 

r)(A;2) = V'vre-K'=^+^'?)'^ (28) 

and the other quantities deflned in Eqs. ([I5]-[18]). Using Eq. (127|) to eliminate (5p/po from 
Eq. (1261) . we arrive at 

ik ■ 5v H — %- (vo ■ 5v) Vq + Vyj5\ — ik [c^ + pu) ^ = —ikfj. (29) 

where, again, = vq — Vp. Since there is no disturbance along the z-direction, we set 
k = k^'k + kyj. Introducing an orthonormal basis {ai, a2, £13} such that 

k . . . . 

ai = k = -, a2 • z = 0, ag = z, (30) 

we project Eq. (!29l) onto these three directions. Deflning 

6v = 6viSLi + 6v2a2 + 6v3a.3, (31) 

Eq. ([29D yields 

2 

ikvr cosi^ip — (p)^^! + z^'^o cos(0) {vq COs(0)(5f i — f sm{(f))Sv2) 



Sv 

+ Vyj^vi - ik{cl + pqu) — — = -ikfi (32) 

Vr COii[ip — (p) 

2 

ikvr cos{ip — (f))6v2 — -^-^Vq sin(0) {vq cos(0)5f 1 — vq sm{(j))6v2) + Vw^V2 = (33) 
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and Svs = 0, where we have changed coordinates from {kx,ky) to {k,(f)). Here is the 
angular orientation of k, measured from the x axis, and ip is the angle between and the 
X axis. We can obtain general results for (5v(x, and 6p{x,y) by solving Eqs. (!32|) . (!33ll . 
and fl27|) for ^Ci, 6v2, and 5p and then inverse Fourier transforming the result. However, for 
simplicity in what follows, we again take r — 0, clamping all of the flocking individuals to 
the same speed. Equations fl32l) . fl33l) . and fl27j) then yield 
6p{k,4>) /c^f/sin^(0) 



Po /i^^t"^ cos^(-?/^ — 0) — /c^c^ sin^(0) — iku^iVr cos^ip 

with 



(34) 



m = c'.+Pom = ci + ^^^. (35) 

By inverse Fourier transforming, we obtain the density perturbation at any point (r, 6) in 
the flock, 

g*fcrcos(<A-9)^^^^^^ (36) 



Po (27r)2 7o 7o po 
where 5p(A;, 0)/po is defined in Eq. f l34l) . In the next section, we explore Eq. fl36l) analytically 
in the case of a static obstacle. After that, we evaluate Eq. (!36|) numerically for both a static 
and a moving obstacle/predator. 



V. ANALYTICAL RESULTS FOR A STATIC OBSTACLE 

To evaluate the integral, we consider the following illustrative case. We take v,. = vq, 

which corresponds to a stationary obstacle or predator. This implies that ip = 0. Also, 

assume that the parameters are such that for most values of k, the quantities c and can 
be approximated by their large k limits. We have 

c~ lim c{k) = Cs, (37) 

fc— >oo 

1- /T\ 4:TTWq 

^ hm Uui{k) = — -— . (38) 

The range over which this approximation is good will be investigated in Sec. IVI[ With these 
approximations we can integrate fl36|) to obtain the density perturbation 6p/po- The full 
analysis is done in the Appendix. Fig. [3] displays the density perturbation for a particular 
choice of the parameters. Note that the main feature is a wedge formed from the information 

11 



of the obstacle propagating through the flock. We flnd, in the Appendix, that density will 
be large when y is near ij/o, where 



1 / /2 
— \ X e 

7 V 27 



(39) 



where we have introduced the quantities e = v^jlc and 7 = fo/c. The first term in Eq. 
( l39l) . x/7, corresponds to the wedge condition, Eq. ( l23l) . 

Figure m shows 5p{x,y)/ po versus y for several fixed values of x. Numerical data (com- 
puted in Sec. IVll) are plotted as open circles, and the theory obtained in the Appendix is 
plotted as a solid curve. They agree well. Further approximations (see the Appendix) result 
in analytic expressions for the height, 7i, and width, W, of these profiles (illustrated in Fig. 
m^a)). The height, Ti {6p/po at maximum), and width, W (distance between maximum and 
minimum), are 

?7o7 /TT 



n{x) 



c2(l + 72) V 2 

2 yrry 



exp 



1 e^P - Ax'je' 
~2 ^ 



(40) 
(41) 



V2 7 

From these expressions we see that the width is predicted to be insensitive to many of the 
parameters of our problem except I and 7. For example, the width does not increase far from 
the source of the disturbance (i.e., W in Eq. (HTj) does not depend on x). A main feature 




FIG. 3: Plot of the density fluctuations. Parameters used for the figure are: Cs = 15, 770 = Ij 
e = 0.838, 7 = 2, and / = 0.1. The x axis is horizontal, the y axis is vertical. 
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FIG. 4: (color online) Plot of the density fluctuation for constant values of x for the same parameters 
in Fig.[3l The red circles are numerical values (computed via Sec. IVII methods) and the solid curve 
is the theory, a) x = 2. The definitions for height and width are displayed on the plot, b) shows 
X = 1.25, and c) x = 0.5. 
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of Eq. (HOj) is its prediction of the exponential decay of the height of the disturbance with 
increasing x. In the next section, we compare numerical simulations with these predictions. 



VI. NUMERICAL RESULTS 

tan 
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FIG. 5: (color online) Plot showing the wedge angle 9 vs. the quantity for = 15. This shows 
the agreement of the numerical data for wq = 0.001 (red circles), wq = 0.01 (green boxes), and 
wq = 0.1 (blue crosses). The solid line is the theoretical prediction for the wedge angle from Eq. 
([21]). The other parameters for the plot are uq = 0.1, Kp = 0.1, and k^, = 0.5, / = 0.1, rjQ = 1, and 
po = 0.8. 



A. The Static Obstacle 

In order to evaluate the integral in Eq. ( |36l) numerically, we need to do a two-dimensional 
infinite integral at each point in physical space. To do this, we express the kernel of the 
inverse Fourier transform as a sum of Bessel functions using the Jacobi- Anger expansion (see. 



for example, [25[)- This allows us to do one of the iterated integrals via contour integration. 
We then obtain an infinite sum of single integrals at each real space point that we evaluate 
numerically. An example of the density fiuctuation evaluated using these methods looks 
very similar to Fig. [31 Similar to the analytic result in the Appendix, the numerical density 
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fluctuation shows a prominent wedge emanating from near tlie origin. Tlie correspondence 
to tlie analytic work is excellent and can be seen in Fig. |H We now compare the numerical 
results to theoretical predictions for 9±^ 7i, and W. 

Using the relation in Eq. fl23l) . we can test the above results to determine the accuracy of 
the numerical fit to the wedge angle predicted earlier. Visually determining the angle from 
the numerical output gives a well defined wedge angle to about 0.5° accuracy. The tangent 
of this angle, Eq. (l23l) . can then be compared with the quantity Cs/vq. Figure [5] shows that 
this comparison yields very good agreement for the static case. The parameters used in Fig. 
|5]are wq = 0.1 (blue crosses), wq = 0.01 (green boxes), and wq = 0.001 (red circles). It is 
seen that changing wq leaves the wedge angle essentially unchanged. 

Figure [6] shows comparisons between results for W from the numerical simulations (col- 
ored markers) with the predictions given in Eq. fj4T]) (solid curves). Figure [6]^a) and shows 
that, as predicted by the theory, W/l is insensitive to the values of Kp, n^, uq, and rjo. The 
only important dependence of the width was on the parameter 7 as seen in Fig. Mjo). Here 
we see that the wedge width approaches a constant value for large 7. Figure [71^a-d) show 
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FIG. 6: (color online) Graphs showing the dependence of the width on various parameters of the 
model. The solid curves is the expression in Eq. (j41l) . The colored markers are numerical values 
obtained using the processes described in Sec. lVIi a) W/l vs. Kp (squares), (circles), uq (crosses), 
and r/o (triangles), b) W/l vs. 7. If a parameter is not varied then it has the value: pQ = 0.8, 
Cs = 15, 7 = 2, r/o = 1, no = 0.001, wq = 0.5, = 0.5, Kp = 0.1, / = 0.1, and x = 2. 
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the dependence of the height, 7i, on /t^,, wo, rjo, and 7, respectively. In these figures, as well 
as in Fig. [6l if a parameter is not varied, then it has the value: po = 0.8, Cg = 15, 7 = 2, 
rjo = 1, uq = 0.001, Wo = 0.5, = 0.5, Hp = 0.1, / = 0.1, and x = 2. As predicted by Eq. 
fj40|) . H is linear in rjo. Figure [7](e) shows that the height is insensitive to both Hp and uq. 
The theoretical prediction for the dependence of the height on position, x, is verified in Fig. 
WIS)- Figure [Tl^g) shows that there is agreement with Eq. PUI) for / ^ 0.15, but breaks down 
at small / since the width is predicted to go to zero in that case. Additionally, the expansion 
in the Appendix used to obtain Eq. flA9|) implies that our approximations are expected to 
become invalid at large e/7 = 2ttwo/ h^vq. For example, at very low k^, the width starts 
deviating from the prediction (Fig. [6](a), circles). 

B. The Moving Obstacle 

We numerically evaluated Eq. fl5B]) using the same method as Sec. IVI A[ but with nonzero 
predator angle, ip. The results of the comparison between the theoretical prediction of the 
wedge angles and the numerics can be seen in Fig. [Si The theory, Eq. ([23]), predicts the 
wedge angles as a function of predator angle, ip. The figure shows the correspondence to 
the numerical data for two angles, ip = 7r/3 and ip = tt/6, versus various values of Cs/vq. 
The agreement is excellent, and, similar to the static case, the wedge angles are insensitive 
to parameters such as the nonlinear viscosity parameters, Wq and k^j- 

VII. CONCLUSIONS 

In this paper we have explored the response of a fiock to static and moving obstacles. The 
obstacle is introduced into a fiock and the fiuctuations about an equilibrium are analyzed. 
We find that with both the static and the moving obstacles, the fiock produces a prominent 
wedge where the information is propagating away from the disturbance, as shown by Fig. 
[21 The wedge angles can be predicted using a simple geometric construction. The infor- 
mation/disturbance propagates asymmetrically (unless tp = 0), with two angles, 6^ and 6'_, 
given by Eq. (^^. We tested this analytically as well as numerically, and the result is found 
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to agree well with the theoretical prediction. The wedge angles are insensitive to most phys- 
ical parameters, most notably the velocity viscosity term, W, and, unlike the well-known 
Mach cone in acoustics, there is no threshold speed for existence. Specifically, the angles 
only depend on the speed of sound in the fiock, c^, the speed of the fiock, vq, the relative 
speed of the obstacle to the fiock, Vr, and the angle between them ip. Heights and widths 
of the Mach cones for ip = are given by the analytic expression in Eq. fHOl) and Eq. fHTj) . 
Numerical results are in good agreement with these expressions. It is also noteworthy that 
the wedge width, defined in Fig. 111(a), is insensitive to many parameters in the model as can 
be seen in Fig. |6](a). 

Future work should include the dynamics of an obstacle hitting a fiock, extension to 
r 7^ 0, and a physical explanation of the wedge shape and offset from the origin. Finally, 
the extension to a fully nonlinear treatment of the obstacle is of interest. 

This work was financially supported by ONR grant N00014-07-1-0734. 
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APPENDIX A: DERIVATION OF 5p/pQ, THE HEIGHT, AND THE WIDTH OF 
THE DISTURBANCE 



We derive the density perturbation Sp/po, the height, H, and the width, W, via a direct 
computation of the integral Eq. fl36l) . To evaluate the integral, we consider the following 
special case. First, we take v,, = vq. This implies that ip = 0. Also, as in the main body of 
the paper, 

c ^ lim c{k) = Cs, (Al) 

fc— >oo 

z/^ ^ lim iy^{k) = (A2) 

and we define e = ^ and 7 = ^. With these approximations we can write the integral (|36|) 
in rectangular coordinates, 

- ' e^fc^^e 4^=. / , , dkydK. (A3) 



Let us do the ky integral first. We shift the path of integration up in the complex ky plane 
to lm{ky) = so as to go through the saddle point in the complex plane giving. 



where the integral is over real u, we have factored the denominator, and we define 



, z2e i2y 

ui = -jk,\ll- — -—, (A5) 



^2 = 7^x^1-;^-^. (A6) 



Using 

/ 7^ ^ -du = (a - iuiYw{bui) + {a - tU2rw{-bu2) \ + ^ (AT) 

J [U-Ui)[u-U2) U1-U2 ' b 

along the contour given in Fig. [9l we can explicitly evaluate the u integral in terms of the 



complex error function w, given by (see 25l |) 



00 —t 



2 



w{z) = - dt = e-''erfc(-2z) if lm{z) > (A8) 

TT ./-oo z - t 
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and defined for tlie negative imaginary z by analytic continuation. 

Expanding Eq. (lASP and Eq. flA6|) as Taylor series in [ej^k^^ we obtain from Eq. flA4 



Po 

where 



= C±(u;)erfc(F±(w)), 



and 



I/O 

5 



1 / /2 

— X e 

7 V 27 



^0^ 



exp 



- — (x 

7 



47 



27 



7 



zW 



" f(l + 7^) 
7(3; - l^e) 



± 



ry 



f(l + 7^) 



f^dw 
(A9) 



(AlO) 

(All) 
(A12) 

(A13) 

(A14) 

(A15) 



VTT^ /(1 + y) /(I + 7^) 
We now approximate the integral over w using, for example, integration formula 25.4.46 



on pg. 890 of 



25|, to obtain an analytic expression for Sp/po. All non-numerical plots and 



images are composed via this method (using n = 10). We can further use 



erf(M + iv) ~ erf('u) 



27ru 



[(1 — cos(2'ut>)) + i sin(2-ut>)] ^ erf(M) 



(A16) 



from pg. 299, of the same text, to approximate the integrand. For large enough x, if we 
change variables and shift the origin in the y direction to the center of the wedge, yo, we see 
that the real part of the argument of the error function is 

where y = yo + y. Since for modest values of x this is typically far from zero, the error 
function is approximately constant and equal to 2. This gives (near the center of the wedge 
for fixed x) 



A{x)e 



' f+{w)dw ~ A{x)e 



47^v/7r 

/2(l + 72) 



(A18) 
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where we have integrated a Gaussian, and neglected the imaginary part, since the final 
integral must be real. Thus, the height, 7i {6p/po at maximum), and width, W (distance 
between maximum and minimum) defined in Fig. HJ are given in the main body of the paper 
(Eq. (go]) and Eq. (gl])) as 



n{x) 
w 



Vol 

c2(l + 72 

2 



1 + 7^ 

7 



2exp 



1 e2/2 
-2 + — 



4^2 
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FIG. 7: (color online) Various graphs showing the dependence of the height of the wedge on various 

parameters of the model. The solid curves is the expressions in Eq. (j40p . The colored markers are 

numerical values obtained using the processes in Sec. IVII a) 7i vs. Kw b) Ti vs. wq. c) 7i vs. r/g. 

d) TC vs. 7. e) 7i vs. Kp (boxes), and uq (circles), f) Ti vs. x. g) Ti vs. I. If a parameter is not 

22 

varied then it has the value: po = 0.8, Cs = 15, 7 = 2, tjq = 1, uq = 0.001, wq = 0.5, Kw = 0.5, 




FIG. 8: (color online) Plot showing the wedge angles 6± vs the quantity This shows the 
agreement of the numerical data for two choices of predator angle tp. The solid lines are the 
prediction of Eq. (j24p for ip = vr/S, whereas the dashed lines are for -0 = 7r/6. The lines with the 
positive slope correspond to 0+ (blue crosses), and the lines with the negative slope correspond to 
0_ (red circles). The lines are the theoretical prediction for the wedge angles. The other parameters 
for these plots were wq = 0.001, uq = 0.1, Kp = 0.1, and = 0.1. 
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FIG. 9: Contour in complex plane for the integral in Eq. (|A7p . This contour is forced by causality 
and from y ^ —y symmetry. 
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